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FORTRAN PROGRAM FOR CALCULATING AXIAL TURBOMACHINERY 

BLADE COORDINATES 
by Theodore Katsanis 
Lewis Research Center 

SUMMARY 

A FORTRAN IV computer program has been written to calculate blade coordinates 
with respect to the true blade chord. The required input is the axial blade chord, blade 
stagger, leading- and trailing- edge radii, angles of tangency on leading- and tr ailing- 
edge radii, and a few intermediate spline points. This input is identical to the geomet- 
rical input required for blade-to-blade aerodynamic analysis programs previously pub- 
lished by NASA (TN D-5427, TM X-1764, and TN D-5044). 

INTRODUCTION 

There are several NASA computer programs for calculating velocities on a blade-to- 
blade surface between blades (refs. 1 to 3). These programs are easy to use because the 
blades can be described very simply. The required geometrical input consists of the 
axial blade chord and stagger, leading- and trailing-edge radii, angles of tangency on the 
leading- and trailing-edge radii, and a few intermediate points which are fitted with a 
spline curve. This required geometrical input results in a precisely defined blade sur- 
face. After a satisfactory blade surface velocity distribution is obtained, it is often de- 
sired to calculate a large number of offset coordinates with respect to the true blade 
chord. The true blade chord is tangent to the lower surface of the blade. Since the blade 
shape is specified by mathematical equations, these coordinates may be calculated in a 
straightforward manner. However, this is a tedious and time consuming hand calcula- 
tion. It is the purpose of the program TFORM to perform these calculations. 

The FORTRAN IV program TFORM is presented herein with a complete description 
of the input required and the output obtained. The input and output for an example case 
are also given. The geometrical input is just a part of that required for the programs 
TSONIC, TURBLE, or TANDEM (refs. 1 to 3). 



SYMBOLS 


r radius from axis of rotation 
w linear coordinate in tangential direction 
w- coordinate of blade surface 
Wq w-coordinate of (x, y) origin 

Wj w-coordinate of (x^ 0) 

x coordinate tangent to blade lower surface 

11. 

x. x- coordinate at 1 increment from blade leading edge 

y coordinate normal to x-axis 

ii. 

y. y- coordinate at i 1 n increment from blade leading edge 

y^ . y for lower blade surface 

y u ^ y for upper blade surface 

z axial distance from blade leading edge 
z- coordinate of blade surface 
Zq z-coordinate of (x, y) origin 

Zj z-coordinate of (x., 0) 

0 angular coordinate about axis of rotation, radians 
c p blade angle from axial direction, deg 


TRANSFORMATION PROCEDURE 

The basic transformation consists of a rotation and translation. The input coordi- 
nates are given as (z, 0) coordinates where z is the axial direction and 9 is the angular 
coordinate in radians about the axis of rotation. The linear coordinate in the 0-direction 
is equal to r0 = w. The origin in the w-z plane is at the leading edge of the blade, as 
shown in figure 1. The entire curve for each surface is specified mathematically by the 
leading- and trailing- edge radii and by a spline curve in between. The output coordinates 
are given as (x, y) coordinates with the x-axis tangent to the blade lower surface and the 
y-axis tangent to the blade leading edge, as shown in figure 2. 

The first step in the program is to determine the angle cp, the true chord, and (z q ,Wq) 
which are the (z,w) coordinates of the (x,y) origin. These constants specify the amount of 
translation and rotation and are calculated by equations (Al) to (A18) in the appendix. 
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w = r6 



y 



The next step after the translation and rotation constants have been calculated is to 
calculate the y- coordinates for each blade surface corresponding to each increment in 
the x-direction. This is done by finding the intersection of the line x = x^ and the curve 
w = w(z) (see fig. 2). hi w-z coordinates, the line x = x. is 

z - z. 

w = W* (1) 

tan cp 
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where 


Z 1 ~ z 0 + x i C0S ^ > 

w* = w n + x. sin cp 
1 U X > 


( 2 ) 


The blade surfaces are described mathematically by piecewise functions; that is, the 
leading- and trailing- edge segments are given by the equation of a circle and the rest of 
the blade by a spline curve which is a piecewise cubic polynomial (ref. 4). We can de- 
note this by 


w = w b (z) (3) 

For any z then, w is determined as indicated by equation (3). Equations (1) and (3) can 
be solved simultaneously to determine (z b , w b ) where the line intersects the blade. The 
numerical procedure for solving equations (1) and (3) is described in the appendix. Then 
y. is calculated by 

y t -f - z b ) 2 + (w x - w/ (4) 

DESCRIPTION OF INPUT AND OUTPUT 

The computer program requires as input a geometrical description in (z, 0 ) coordi- 
nates of the two blade surfaces, the radius r, a scale factor if desired, and the desired 
x- increment for the output coordinates. Output from the program includes x- and 
y- coordinates for the upper and lower surfaces (see fig. 2). 


Input 

Figure 3 shows the input variables as they are punched on the data cards. The first 
input card is for a title, which will serve for problem identification. The remaining 
cards are for input variables. All variables are real (decimal point must be punched) in 
a 10- column field. It should be noted that the input corresponds very closely to the blade 
geometry input for the NASA blade-to-blade analysis programs of references 1 to 3. 
Further explanation of the input variables is given in the Instructions for Preparing Input 
section. 
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1 10|11 20(21 30|3l 40(41 50(51 60(61 70|71 80 

TITLE 


CHORD 

STGR 

RMI 

SCALE 

DELX 

I^yS.Z I * 






RI1 

ROl 

BET 11 

BETOl 

SPLNOl 






MSP1 ARRAY I 

















THSP1 ARRAY i 


















RI2 

R02 

BETI2 

BET02 

SPLN02 


jp— 






wmimm. 

MSP2; 

\RRAY | 

















THSP2 ARRAY 1 


















Figure 3. - Input form. 


The input variables are as follows: 


CHORD 

STGR 

RMI 

SCALE 

DELX 


RI1, RI2 

ROl, R02 

BETI1, 

BETI2 

BETOl, 

BET02 


Overall length of blade in the z-direction, see fig. 4 

Angular 0-coordinate for center of trailing-edge circle of blade with respect 
to center of leading- edge circle, radians, see fig. 4 

Radius of blade section from the axis of rotation (If RMI = 1, then all 
0- coordinates are the actual linear dimension w. ) 

Ratio of output dimensions to input dimensions (For example, if input is in 
feet and output is desired in inches, SCALE = 12 should be used.) 

Spacing of output coordinates in the x-direction, see fig. 5 (DELX should 
be chosen to be at least CHORD*SCALE/100. DELX must be given in the 
output units; i. e. , if input is in feet and output is in inches (SCALE = 12), 
then DELX is in inches. ) 

Leading-edge radii of the two blade surfaces, see fig. 4 

Trailing-edge radii of the two blade surfaces, see fig. 4 

Angles (with respect to z-direction) at tangent points of leading-edge radii 
with the two blade surfaces, deg, see fig. 4 (These must be true angles 
in degrees. ) 

Angles (with respect to z-direction) at tangent points of trailing-edge radii 
with the two blade surfaces 
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SPLNOl, 

SPLN02 

MSP1, 
MSP 2 


THSP1, 

THSP2 


Number of blade spline points given for each surface as input, maximum 
of 50 (These include the first and last points (dummies) that are tangent 
to the leading- and trailing- edge radii (fig. 4).) 

Arrays of z-coordinates of spline points on the two blade surfaces, meas- 
ured from blade leading edges, see fig. 4 (The first and last points in 
each of these arrays must be left blank, since these values are calculated 
by the program. If the last point is on a new card, a blank card must be 
used. ) 

Arrays of 0-coordinates of spline points corresponding to MSP1 and MSP2, 
radians, see fig. 4 (Blanks must be used in positions corresponding to 
those in MSP1 and MSP2. ) 


Instructions for Preparing Input 

Units of measurement . - Two units are used: one for linear measurements and one 
for angles. Any unit may be used for linear measurement. If a different unit is desired 
for output, this may be accomplished by the use of a scale factor in SCALE. If SCALE 
= 1, the output units are the same as the input units. The angular measurement 0 must 
be given in radians. However, if RMI = 1 is specified, the 0-coordinate can be given 
as a true linear measurement. 

Blade geometry . - The upper and lower surfaces of the blade are each defined by 
specifying three things: leading- and trailing- edge radii, angles at which these radii are 
tangent to the blade surfaces, and z- and 0-coordinates of several points along each sur- 
face. These angles and coordinates are used to define a cubic spline curve fit (ref. 4) to 
the surface. The standard sign convention is used for angles, as indicated in figure 4. 

The blade must be oriented with a concave lower surface. 

A cubic spline curve is a piecewise cubic polynomial which expresses mathematic- 
ally the shape taken by an idealized spline passing through the given points. Reference 4 
describes a method for determining the equation of the spline curve. When this method 
is used, only a few points are required to specify most blade shapes accurately, usually 
no more than five or six, in addition to the two end points. As a guide, enough points 
should be specified so that a physical spline passing through these points would accurately 
follow the blade shape. This means that the spline points should be closer where there is 
large curvature and farther apart where there is small curvature. 

The coordinates for either surface of the blade are given with respect to the leading 
edge, with the leading edge of the blade being defined as the furthest point upstream. 

Format for input data . - All input variables are real numbers (punch decimal point) 
in a 10- column field. 
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Figure 4. - Geometric input variables. Angles BETH, BETI2, BETOl, and BET02 must be given as true angle 
in degrees, notangle as measured in z-0 plane. 


TABLE I. - INPUT FOR SAMPLE PROBLEM 


2ND ROTOR 

HUB SECTION 




CHORD 

STGR 

RMi 

SCALE 

DEL X 

0® 7966700E-01 

- 0® 36150QQE-0I 

0« 259 3300 

I 2 *000000 

0 * 5Q0OOOOE— 01 


BLADE SURFACE 

1 — UPPER SURFACE 



Rii 

RG1 

3E Ti 1 

BETOl 

SPLN01 

0 ® 31250QQE— 02 

0* 8330000E-03 

41 « 000000 

-46 .300000 

3*0000000 

MSP1 AR< AY 





0 

0« 5883 OQOE— 01 

-0 



THSP1 ARRAY 





0 

0® 345400 GE— 01 

-0 




BLADE SURFACE 2 

— LOWER SURFACE 


R I 2 

R02 

8ETI2 

BET02 

0*31250006-02 

0* 8330000E— 03 

29« 000000 

-35 * 500000 

MSP 2 ARRAY 




-0 

0*58830006-01 

-0 


THSP2 ARRAY 




-0 

0* 4820 000E-02 

-0 



SPLMQ2 

3*0000000 
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TABLE n. - OUTPUT FOR SAMPLE PROBLEM 


BLADE DATA AT INPUT SPLI NE POINTS 


2 

0.107486- C2 
0 • 58830E— 0 1 
0.79436E-C1 


BLADE SURFACE 1 


THETA 

0.90945E-02 

0.34540E-01 

-0.33931E-01 


DERIVATIVE 
3.35205 
-2.45348 
-4. 03517 


2ND DERIV. 
-101.424 
-99.6146 
-53.9014 


l 

0 #464006— 02 
0 » 58830E- 01 
G * 78 35 0E — C 1 


BLADE SURFACE 2 


THETA 

-0*1 0539E- 01 
0.482 00E-02 
-0.38765E-01 


DERI VATI VE 
2.13747 
-1.62164 
-2. 75052 


2ND DERIV. 
-66.5432 
-72.1948 
-43.4682 


NO# OF POINTS = 21 PHI = -5*3371 DEGREES 


X 

Y LOWER 

Y UPPER 

0 

0. 3 75 00E— 01 

0. 375 00E— 0 1 

0# 5000 0E- Cl 

0.21445E-02 

0.10269 

0.10000 

0.3332 5E-01 

0.14606 

0.15000 

0.61872E-01 

0.18273 

0.20000 

0. 86266E-01 

0.21294 

0.25000 

0.1065 7 

0.23691 

0.30000 

0. 122 85 

0.25485 

0.35000 

0. 13516 

0.26696 

0.40000 

0. 14356 

0.27343 

0.45000 

0. 14812 

0.27443 

0.50000 

0. 14888 

0.27014 

0.55000 

0. 14591 

0.26071 

0.60000 

0.1392 6 

0.24630 

0*65000 

0.12899 

0.22705 

0.7000 0 

0.11513 

0.20310 

0.75000 

0.978Q8E— 01 

0.17468 

0.80000 

0. 77298E-01 

0.14224 

0.85000 

0. 53936E-01 

0.10626 

0.90000 

G.28038E— 01 

0.67177E-01 

0.95000 

0.32655E-03 

0. 25395 E— 01 

0.96253 

0* 99960E — 02 

0. 99960E— 02 


Output 


Sample output is given in table II for the example blade given in table I. The first 
output gives additional computed blade data at the input spline points. This includes the 
z- and 0-coordinates at the points where the spline curves are tangent to the leading- and 
trailing- edge radii. Also, the first and second derivatives are given at each spline point. 
Of particular interest are the second derivatives. Any error in blade geometry input will 
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DELX 



Figure 5. - Output coordinates. 


usually result in wild values for some of these second derivatives. 

The next output gives the transformed blade coordinates. The first line of output 
gives the number of x-coordinates and the orientation angle cp, as shown in figure 5. 
This is followed by a tabulation of the x- and y-coordinates for the upper and lower sur- 
faces (see fig. 5). 


Error Conditions 

The error message is given first for each error condition: 

(1) BETI2 MUST BE GREATER THAN PHI AND BET02 MUST BE LESS THAN PHI 
TO HAVE X AXES TANGENT TO LOWER BLADE SURFACE 

It is assumed in the program that the x-axis is tangent to the leading- and tr ailing- 
edge radii. If either BETI2 is less than cp or BET02 is greater than cp this tangent line 
will not actually be tangent to the lower blade surface, and part of the lower surface will 
be below the x-axis. Normal calculations will still be made, but there will be negative 
values for Y LOWER. 

(2) PART OF BLADE HAS NEGATIVE X VALUES 

This message is printed if part of the blade would extend to the left of y-axis. This 
can happen if BETH is greater than <p + 90° or if BETI2 is less than cp - 90°. No fur- 
ther calculations are made and the program will proceed to the next case. 

(3) LOWER BLADE SURFACE IS NOT ENTIRELY CONCAVE 

This message is printed if some part of the blade lies below the x-axis. Normal 
calculations will still be made, including negative values for Y UPPER or Y LOWER. 

(4) Z COORDINATE IS NOT WITHIN BLADE 

This message is printed by subroutine BLCD if the z-coordinate given this sub- 
routine as input is not within the bounds of blade surface. The value of z and the blade 
surface number are also printed when this happens. This message should only occur if 
there is an error in the input data. 
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(5) ROOT HAS FAILED TO OBTAIN A VALID ROOT 

This message is printed by subroutine ROOT if a root cannot be located, or if the 
accuracy of the root is not satisfactory. The user should thoroughly check the input data. 


PROGRAM PROCEDURE 

The main program is TFORM. There are 4 subroutines: FUNCT, ROOT, BLCD, 
and SPLN22. The calling relation of all the subroutines is shown in figure 6. 


TFORM 

’ 

L 



| root] 




r 



FUNCT | 


E 

| BLCD | 

□ 


‘ SPLN22 


Figure 6. - Calling relation 
of subroutines. 


TFORM reads and prints out all the input data. Then the transformation constants 
cp, z Q , and w Q are calculated as described in the appendix. Next the x and y arrays 
are calculated. The method for calculating y for a given x value is described in the 
appendix. The root finding procedure required by this method is accomplished by sub- 
routine ROOT. 

Subroutine FUNCT calculates f(z) in equation (A19) for either the upper or lower 
surface. 

Subroutines ROOT, BLCD, and SPLN22 are the same as described in references 1 
to 3. Subroutine ROOT was changed in reference 1 from the coding used in references 2 
and 3. This was to adopt the more foolproof method of locating roots by the bisection 
method. Subroutine BLCD calculates the 9 blade coordinates when given a z-coordi- 
nate. Subroutine SPLN22 calculates the spline curve for the blade surfaces. 
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FORTRAN Variables in TFORM and FUNCT 


A 

B 

BETI 

BETO 

C 

CHORD 
CNVX 
CP HI 
D 

DELX 

E 

F 

FZ 

G 

H 

I 

INDEX 

ISURF 

J 

MSP 

NOPT 

NSP 

NSPI 

P 

PHI 

PHICC 

PHICOR 

PHIDEG 


a, fig. 7 

b, fig. 7 

array, BETI1 or BETI2, see input 
array, BETOl or BET02, see input 

c, fig. 7 
see input 

if CNVX > 1, either the lower or the upper surface has negative y- coordinates 
cos cp 

d, fig. 7 
see input 

e, fig. 8 

f, figs. 8 and 9 
f(z), eq. (A 19) 

g, figs. 8 and 9 

h, figs. 8 and 9 
temporary index 

used as both a switch and subscript in calculating y blade coordinates 

index indicating blade surface number 

index for DO loop 

input arrays MSP1 or MSP2 

number of points in x and y output arrays 

number of spline points 

array of number of spline points 

p, figs. 8 and 9 

cp, fig. 7 


<P C _ C , fig- 7 

<Pcorr’ fi S- 7 
cp, deg 
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PI 


7T 


RI array, RI1 or RI2, see input 

RMI see input 

RO array, ROl or R02, see input 

SCALE see input 

SPHI sin cp 

SPLNO either SPLNOl or SPLN02, see input 

SRW integer code variable that causes either ROOT (if SRW = 21) or SPLN22 (if 

SRW = 18) to write out data useful for debugging 

STGR see input 

TCHORD true chord, fig. 2 

THETA d 

THSP input arrays THSP1 or THSP2 

TOLERW permissible tolerance in value of w for a given value of x, TOLERW 
= CHORD x 10“ 4 

TP HI tan <p 

WO Wq, fig. 2 

WB w b , fig. 2 

W1 w p fig. 2 

X array of output values of x 

Y array of output values of y^ and y u 

Z0 z Q , fig. 2 

Z1 z p fig. 2 

ZB z b , fig. 2 

ZERO zero value variable 

ZL, ZT if Z1 is less than ZL or greater than ZT, the blade surface in the x, y- 
coordinates is the opposite of the one in the w, z-coordinates 
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PROGRAM LISTING 


c 

c 

c 


c 

c 

c 


COMMON SR W t INI T(2) *TPH I *wB ,Zi ,W1 tCHORQ* STGR , RMI * RI ( 2) , RQC2 )* 
1 6ETUZ) ,BhTG<2) ,NSPI (2) ,MSP<50,2) ,THSP(50,2) 

DIMENSION X(10I),Y(iQl,2) 

REAL MSP 

EXTERNAL FUNC T1 ,F UNC T2 
1 CONTINUE 
1NITC 1) = 0 
IN I T{ 2) = 0 

READ AND PRIHJ ALL INPUT DATA 


WRITE 
READ! 
WRI TE 
WRI TE 
READ 
WRI TE 
DO 10 
IF (J 
IF (J 
WRITE 
READ 
WRI TE 
NSP K 
NSP = 
WRITE 
READ 
WRITE 
WRI TE 
READ 
10 WRI TE 


{ 6, IOOOI 
5*1100) 

{ 6 * 1 100 ) 

{ 6*1110) 

( 5*1030 ) CHORD * STGR * RMI * SCALE* DEL X 
16,1040) CHORD, STGR, RMI , SCALE , DELX 
J=l * 2 

» EU. 1 ) WRI TE (6,1120) 

.EQ. 2 ) WR I TE (6 * 1 130) 

(6,1140) J,J*J,J,J 

( 5,1030) RI( J) ,RO(J) ,B£ T I ( J) * BETO ( J ) 
(6*1040) RI(J),RO(J) *BETI ( J) * BETOC J) 
J) = SPLNO 
NSP I ( J) 

(6,1150) J 

( 5*1030) (MSP (I , Ji ,1=1 , NSP) 

( 6,1040) (MSP (I , J) ,1 =1 ,N$P) 

(6*1160) J 

(5,1030) { THSP(I , J) ,1=1, NSP) 

( 6,1040) ( IHSPCI , J) ,1=1, NSP) 


, SPLNO 
* SPLNO 


CALCULATE TRANSFORMATI ON CONSTANTS 


PI = 3.1415927 
CNV X = 0. 

TOLERW = CHORD /l. E 4 
A = CHQRD-RI(2)-R0<2) 

B = STGR*RMI 
C = SQ RT( A*A+B*B ) 

PHICC = A TAN ( 8/A ) 

PHICOR = ARSINi CR£(2)-R0(2))/C) 

PHI = PHI CC+PHICQR 
PHIDEG = PHI /P I * 1 80* 

IF( BETH2) .LT.PHIDEG.0R.BETQ(2) .GT. PHIDEG) WRITE (6,1165) 
IF(BETl(l)~90..LE. PHIDEG. AND. 90. + BE Ti ( 2 ) . GE* PH I DEG ) GOTO 15 
WRITE(6,1167) 

GO TO 1 
15 CGNTINLE 

SPHl = SIN(PHI) 

CPH I = CO S( PH I ) 

TP H I = TAN { PH I ) 

D = C*C0S(PHIC0R) 

E = ( R I ( 1 )-RI (2) )*CPHI 
F = R IC2)*SPHI 
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G = (RUlJ-E)*CPHl 
H = R 1(2) *CPH I 
P = ( R I ( 1 ) — £ ) * SPH I 
TCHORO = U-E+RI iU+ROiZ) 

IF( PH I <*Lt « 0* ) GG TO 20 
E - ( RG( U-RQ(2) i*CPHl 
G - R I (2)*CPHI 
P = F 

TCHORO = 0-E*RI42URGUJ 
20 ZO = Rl(2)+F-G 
WO = -H-P 

CALCULATE X AND V ARRAYS 
XU) = 0. 

Y 4 1 , 1 ) = R I ( 2 ) * SCALE 

IF ( PH I *LT • 0* J YU, U = Y*4IUUSPHI*fRI (2)-RI(l) )*SCAL£ 
Y( 1, 2) = Y(iil) 

I = 1 

35 I = l + l 

XU ) = X( I-D+DELX 
ZI = ZC+XII) /SCALE*CPHi 
W1 - WQ+XC I)/SCAL£*SPHI 
ZERO = 0. 

ZL = -Wl* TPH I 

IT = CHORD- 1 Wi-RHI*STGR)*TPHI 
OU IOC I S URF = 1 ? 2 
INDEX = I SURF 
A * 0. 

8 = CHORD 
ZB = Z 1 

I F ( PH I *EQ « 0* ) GO TO ( 60,80) ,1 NDEX 
IF ( PH I* i FLOAT 4 I SURF ) - 1 ® 5) « LE® 0« ) GO TO 40 
C PHI NEGATIVE AND UPPER SURFACE OR 

C PHI POSITIVE AND LOWER SURFACE 

A = RI( IN 0£X)*(1®-C PH I ) 

IF{ Zl ® LE * Z T) GG TO ( 5 0 , 70) ? I NDE X 
INDEX = 3-INDEX 
A = CHORD-ROC INQEX)*< l.-CPHI ) 

GO TO (50,70) , INDEX 

C PHI NEGATIVE AND LOWER SURFACE OR 

C PHI POSITIVE AND UPPER SURFACE 

408= CHORD-ROC INOEX)*( l.-CPHI ) 

i F ( Z 1 * GE ® ZL ) GO TO C 5 0 , 70) , I NOE X 
INDEX = 3-INDEX 
B = R IUNDEX)*4 l.-CPHI ) 

GO TO (50,70) , I NDEX 

50 CALL ROOTCAtB , ZERO, FUNC TI , TULERW ,Z 8) 

60 CALL FLNCT1I Zb ,FZ) 

GO TO SO 

70 CALL ROOT l A, B , ZERO , FUNC T2 , TOLE R W ,Z8) 

80 CALL FUNCT2(ZB,FZi 

90 YU, I SURF I = SQRT ( ( ZB— Zi ) + ( WB— WI ) # $2 ) *SCAL E 

IFS wB « GE* Wl ) GO TO 100 
CNVX = CNVX*L» 

Y ( I , I SURF ) = -YU * I SURF) 

100 CONTINUE 

IFC X( I )+DELX.LT.TCHORD*SCALE. AND.LLTUOO) GO TO 35 
1F( CNVX.GT.O® ) WRITE (6,1190) 
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NOPT = 1+1 

X(NOPT) = TCHORO* SCALE 
Y ( NOP T » 1) = RO ( 2) * SCALE 

IFIPHI.GT.O. ) Y (NOPT , l ) = Y (NOP T , 1 ) +SPHI * ( RO ( 1 ) -RO! 2) )*SCALE 
Y ( NOP T » 2) = Y ( NOP T » 1 ) 

PHI = PHI /P 1* 180. 

WRITE! 6,1170) NOPT, PH I 

WRITE! 6,1180) (X! I) »YU ,2) ,Y! I ,1 1 ,1=1 , NOPT) 

GO TO 1 

1000 FORMAT ( lHi) 

1030 FORMAT (8F10.5) 

1040 FORMAT l 1X,8G16.7) 

1100 FORMAT ( 8 OH 

1 ) 

1110 FORMAT i 5X » 5HCH0R0 »12X»4HSTGR»13X»3HRMI , 12 X, 5 HSC AL E , 12X , 4H0ELX i 
1120 FORMAT (39HL BLADE SURFACE 1 — UPPER SURFACE) 

1130 FORMAT (39HL BLADE SURFACE 2 — LOWER SJ RF ACE ) 

1140 FORMAT ! 7X,2HR 1 , 1 1 ,12X,2HRO,1 1 ,12X,4HBETI , 1 1 , 1 IX ,4HBET 0, II ,11X,5HS 
lPLNO.Il) 

1150 FORMAT ! 7X > 3HM SP» ll»2X»5H ARRAY) 

1160 FORMAT ( 7X,4HTHSP»ll»2X»5H ARRAY) 

1165 FORMAT (111HL BET 12 MUST BE GREATER THAN PHI, AND BET02 MUST BE LE 
1SS THAN PHI TO HAVE X AXIS TANGENT TO LOWER BLAD-E SURFACE/ 1H.) 

1167 FORMAT ( 37HL PART OF BLADE HAS NEGATIVE X VALUES) 

1170 FORMAT {18H1 NO. OF POINTS = ,1 4 ,10X,5HPHI =,G12.4,8H DEGREES) 

1180 FDRMAT 146HL X Y LOWER Y UPPER/ 

1 !2X,3( G13. 5,5X) ) ) 

1190 FORMAT (45HL LOWER BLADE SURFACE IS NOT ENTIRELY CONCAVE) 

END 


SUBROUTINE FUNCT 

COMMON SR W, INI T!2) ,TPHI ,WB,Z1 »Wl ,CH ORD , STGR, RM I , RI ( 2) ,R0!2 ), 
1 BET I! 2 ) ,BET0!2) ,NSPI 12) ,MSP!50,2) ,THSP(50,2) 

ENTRY FUNCT1!Z»FZ) 

CALL 8L1(Z, THETA) 

GO TO 10 

ENTRY FUNCT2! Z,FZ) 

CALL BL2!Z, THETA) 

10 WB = THETA*RMI 

IF! TPHI.NE.O. I FZ = WB- W1 + ! Z-Zl ) /T PHI 

RETURN 

END 


SUBROUTINE BLCO 

BLCD CALCULATES BLADE THETA COORDINATE AS A FUNCTION OF M (=Z FOR AXIAL) 

COMMON SR W, IN I T12) , TPHI , WB » Z1 ,W1 ,CH ORD , STGR, RMI , RI < 2) , R0( 2 ) , 

1 BET I! 2), BE TO! 2) ,NSPi 12) , MSP (50, 2) ,THSP!50,2) 

DIMENSION EM ( 50,2) , AAA (50) 
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INTEGER SR W» SURF 
ENTRY BL 1 1 M» THE TA I 
REAL M,MSP»MSPMM»MMMSP 
SUR F= 1 
S IGN = 1 . 

GO TO 10 

ENTRY BL2(M,THETA) 

SUR F = 2 
SIGN=-1. 

10 CONTINUE 

N SP = NSP I { SURF) 

IF ( INITI SURF ) .EQ. 13) GO TO 30 
INI TI SERF )= 13 

INITIAL CALCULATION OF FIRST AND LAST SPLINE POINTS ON BLADE 

AA = 8ETI ( SURF )/57. 295779 
AA = SIN(AA) 

M SP ( 1 * SUR F ) = RI (SURF)*(1.-SIGN*AA) 

BB = SORT ( i.-AA**2) 

THSP ( 1 1 SURF ) = SIGN*BB#RI ( SURF) /RMI 
BETH SURF ) = AA/BB/RMI 
AA = 8ET0(SURF)/57. 295779 
AA = SIN(AA) 

MSP(NSP»S UR F ) = CHORD-ROI SURF)* (i.-»-SIGN*AA) 

BB = SORT ( l.-AA**2) 

THSP (NSP, SURF ) = STGR + SI 3 N*BB*RO(SURF ) /RMI 
BETOI SURF ) = AA/BB/RMI 

CALL SPLN22CMSPI 1 , SURF ) , THSP! 1, SURF) , BETI (SURF) , BETOI SURF) ,NSP, 
1 AAA, EMI 1, SURF)) 

IF (SURF.EQ.l) WRITE(6,1000) 

HR I TEI 6,1010) SURF 

WRITE (6, 1020) (MSPIIA, SURF) ,THSP(IA, SURF), AAAIIA), EMI I A, SURF), 
1 IA = i,N SP ) 

BLADE COORDINATE CALCULATION 
30 KK = 2 

IF IM. GT. MSP! 1, SURF) ) GO TO 50 
AT LEADING EDGE RADIUS 
IFIM.LT.O.) GO TO 90 

THETA = SORT I M* ( 2 . *RI I SURF) — M) ) *SIGN 

IF I THETA .EQ.O. ) GO TO 40 

RMM = R 1 1 SURF )-M 

THETA = THETA /RMI 

RETURN 

40 THETA = 0. 

RETURN 

ALONG SPLINE CURVE 

50 IF (M.LE.MSPIKK, SURF) ) GO TO 60 
IF IKK. GE .NSP) GO TO 70 
KK = KK+1 
GO TO 50 

60 S= MSP IKK » SURF )-MSP(KK-l , SURF) 

EMKM 1= EM (KK-1, SURF) 
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EMK, = EM{KK, SURF) 

MSPMM= MSP!KK, SURF )— M 
MMMSP = M-MSP! KK— 1 t SURF ) 

THKk= THSP (KK, SURF )/S 
THKM1 = THSP(KK-1, SURF) /S 

THE TA = £MKMi*MSPMM**3/6. /S + EMK*MMMSP**3/6./S + ( T HK-EMK*S/6. )* 
1 MMMSP + { THKM1— EMKM1 *S/6. ) *MSPMM 
RETURN 

AT TRAILING EDGE RADIUS 

70 CMM = CHQRD-M 

1FI CMM .LT. -CHORD/I. E5) GO TO 90 
CMM = A MAX 1 1 0. » C MM I 

THE TA= SQRT! C MM* ( 2. *RQ ( SURF I —CMM ) I *SI GN 

IF l THETA. EQ. 0. ) GO TO 80 

RMM = RO( S URF }— CMM 

THETA = S TGR + THETA/RMI 

RETURN 

8 0 THETA = STGR 
RETURN 

ERROR RETURN 

90 WRITE! 6,1030) M,SURF 
THETA = 0. 

RETURN 

1000 FORMAT ( 1H1 , 1 3X , 33HBLADE DATA AT INPUT SPLINE POINTS) 

1010 FORMAT! 1HL , 1 7X , 16HBLADE SURFACE, 14) 

1020 FORMAT !7X , 1HZ ,1 OX ,5HTHE TA , 10X , 10HDER I VAT I VE , 5X, 10H2ND DERIV. / 
1 (4G15.5) ) 

1030 FORMAT (33H Z COORDINATE IS NOT WITHIN 8LADE/4H Z =,G14.6,10X, 

1 6HSURF = * G 14. 6 ) 

END 


SUBROUTINE ROOTIA ,B , Y , F UNCT , TOLERY ,X) 

ROOT FINDS A ROUT FOR (FUNCT MINUS Y) IN THE INTERVAL ( A, B) 

COMMUN SR W 
INTEGER SR W 

IF ISRW.EQ.21) WRITE <6, 1000) A, B,Y, TOLERY 
XI = A 

CALL FUNC T! X1,FX1) 

IF! SRW.EQ.21) WRI TE (6 , 1010) Xl.FXl 
X2 = B 

10 DO 30 1=1,20 
X = (X1+X2I/2. 

CALL F LNC T ! X,F X ) 

IF! SRW.EQ.21) WRI TE!6, 1010) X,FX 
IF! !FX1-Y)*(FX-Y) .GT.O. ) GO TO 20 
X2 = X 
GO TO 30 
20 XI = X 
FX1 = FX 
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30 CONTINUE 

IF!ABS!Y-FX).LT.TOLERY) RETURN 
WRITE! 6, l 020) A ,3 , Y,X,FX 
RETURN 

1000 FORMAT (32H1INPUT ARGUMENTS FOR ROOT — A =G1 3 .5 , 3X , 3HB = ,G13.5, 
1 3X » 3HY = ,G13.5,3X,8HT0LERY = , G13 . 5/ 16X , 1HX , 17X , 2 HFX J 

1010 FORMAT (8X»G16.5»G 18.51 

1020 FORMAT (37HLRQ0T HAS FAILED TO OBTAIN VALID R00T/4H A = »G14.6» 

1 1 0 X » 3H B =,G14.6,10X,3HY =,G14.6,3HX =»G14.6»4HFX = ,G14.6) 

END 


SUBROUTINE SPLN22 ! X, Y , Yt P , YNP, N, SLOPE, EM) 

SPLN22 CALCULATES FIRST AND SECOND DERIVATIVES AT SPLINE POINTS 
END CONDITION - DERIVATIVES SPECIFIED AT END POINTS 

COMMON SR W 

DIMENSION X(N) , YIN) ,EM(N) ,SLOPEIN) 

DIMENSION SB 1 100) tGI 100) 

INTEGER SR W 
SB! 1) = . 5 

F = ( Y!2)-Y! 1) ) /(X(2)-XI1) J-Y1P 
G( 1 ) - F * 3 . / ( XI 2) — XI l ) ) 

NO=N- 1 

IFIN0.LT.2) GO TO 20 
DO 10 1=2, NO 
A = I XU)-X( I - 1 ) ) / 6. 

C = I XI I+il-XI 1 )) /6. 

W = 2.*IA+C)-A*S8( I— 1 ) 

SB ( I ) = C /W 

F = (Y(l+l)-YH))/(X(I + l)-X(I))-IYlI)-Y(I-l))/lX(I)-X(l-iJ) 

10 Gil) = <t— A*Gl 1-1) ) /W 

20 F = YNP-I Y(NI-Y(N-lll/IXtN)-X(N-l)) 

W = ( XIN)-X(N-i) )/6.*!2.~SB!N-l) ) 

EMIN) = l F-I XINI-XIN-l) ) *G (N-l) /6. ) /W 
DO 30 I =2 , N 
K = N + l-I 

30 EMIK) = G !K)-SB!K)*EMlK+l) 

SLOPE! 1) = l X 1 1 ) — XI 2) )/6.*!2.*EMU)+EMl2) )+IY(2)-Y( 1) )/(XI 2 )-X I 1)) 
DO 40 1=2, N 

40 SLOPE! I ) = IXt I )-X( 1-1 ) ) /6.*12.*EM< I ) +EMI I — II ) +1 Y! I )-Y ( 1-1 ) )/ 

1 IXI D-XII-i) ) 

IF! SRw.EQ. 18) WRITE 16,1000) N, I XI I ) , YII ) , SLOPE! I ) t EM I I ) , I = 1,N ) 
RETURN 

1000 FORMAT 12X.15HN0. OF POINTS = , I 3/10 X ,1HX , 19X , 1HY , 19 X, 5 HSLOPE , 1 5X, 
12HEM/1 AG 20. 8) ) 

END 


Lewis Research Center, 

National Aeronautics and Space Administration, 
Cleveland, Ohio, May 20, 1970, 

720-03. 
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APPENDIX - EQUATIONS FOR CALCULATING TRANSLATION 
AND ROTATION CONSTANTS 

The following equations can be obtained by referring to figure 7: 


a - CHORD - r.(2) - r Q (2) 

(Al) 

c =^/a 2 + b 2 

(A2) 

tan <p = ^ 

(A3) 

a 


r,(2) - r 0 (2) 


SU1 ^corr = 

(A4) 

^ ^ c-c + ^corr 

(A5) 


The angle cp is the desired angle between the x-axis and the z-axis (fig. 2). The amount 
of translation (z^ and Wq) is obtained next. There are two sets of equations, depending 
on whether <p is positive or negative. 



Figure 7. - Quantities required to compute <p. 
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When cp is negative or zero, the following equations hold (see figs. 7 and 8): 


d = C GOS <? corr 

(A6) 

e = |r.(l) - r.(2)J cos cp 

(A7) 

f = r^(2) sin cp 

(AS) 

g = jr.(l) - ej cos cp 

(A9) 

h = r .(2) cos (p 

(A10) 

P = [rj(l) - ej sin cp 

(All) 

z 0 = r i^ + f ‘ K 

(A12) 

w Q = -h - p 

(A13) 

True chord = d - e + r^(l) + Tq(2) 

(A 14) 


20 



Figure 9. - Quantities required to compute true chord for positive (p (e is 
at trailing edge). 

On the other hand, when cp is positive, the equations for e, g, p, and the true 
chord change as follows (see fig. 9): 


6 = [V 1 ) “ V 2) ] cos ‘p 

(A15) 

g = r.(2) cos (p 

(A16) 

p = r^(2) sin cp 

(A17) 

True chord = d - e + ^(2) + r Q (l) 

(A 18) 


This completes the determination of the transformation constants. 


CALCULATION OF y-COORDINATE OF BLADE SURFACE 

The problem here is to find the y- coordinate for a given blade surface corresponding 
to a given x-coordinate. This can be done by the simultaneous solution of equations (1) 
and (3) . If w is eliminated, z can be obtained by finding a root for the function 


f(z) = W b (z) - Wj 



tan cp 


(A 19) 
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There are two problems that arise. First, the function w^(z) must be the correct blade 
surface. This is not straightforward, since y u could be on the lower surface in the 
w-z plan'e near the leading edge. The second problem is that for certain values of z 
and w near the leading or trailing edge there may be two solutions to equation (A19). 
Both of these problems are overcome by restricting the interval for z for which the so- 
lution is found. After the proper interval for z has been determined, the proper sur- 
face can be ascertained so that equation (A19) must have a unique root. 

With the proper interval for z and the correct blade surface the unique root for 
equation (A19) is found by the bisection method. That is, the interval is bisected to 
determine z n , then f(z n ) is calculated to determine whether the root is in the right or 
the left interval. This gives a reduced interval. The procedure is repeated until the 
root has been located within the desired accuracy. 
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